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Abstract 

This paper considers a special case of the problem of identifying a 
static scalar signal, depending on the location, using a planar network of 
sensors in a distributed fashion. Motivated by the application to monitor- 
ing wild-fires spreading and pollutants dispersion, we assume the signal 
to be Gaussian in space. Using a network of sensors positioned to form 
a regular hexagonal tessellation, we prove that each node can estimate 
the parameters of the Gaussian from local measurements. Moreover, we 
study the sensitivity of these estimates to additive errors affecting the 
measurements. Finally, we show how a consensus algorithm can be de- 
signed to fuse the local estimates into a shared global estimate, effectively 
compensating the measurement errors. 



1 Introduction 

In this paper we consider a distributed parametric estimation problem: we 
want to identify a static scalar field, depending on the location, using a planar 
network of sensors in a distributed fashion. We assume the field to be Gaussian 
in space: the goal of the sensor network is to identify the four parameters 
describing the Gaussian function. This assumption can be seen as a simplified 
model of pollutants dispersion or of wild-fires spreading, which the network is 
required to monitor. The objective of this paper is to design a network able to 
achieve such identification goal, using noisy punctual measurements and local 
communication. We are going to design the spacial deployment of the sensors, 
their communication topology, and the algorithms both to infer the parameters 
from the measurements, and to aggregate local information into a common 
estimate. 
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Sensor networks and distributed estimation are a vast and rapidly growing 
field: among the available literature, we refer the reader to [15l [H [6l [5] for some 
recent research works, and to [HI [HI [7| for some broader discussions. 

Among such wide literature, three recent works [HI [IHl [H] are especially 
relevant to us, as they consider the problem of estimating a scalar planar field 
through point-wise noisy measurements, taken by coordinated moving sensors. 
Let us discuss their contributions. First, the authors of [9] consider the problem 
of finding the leading coefficients of a suitable series expansion of the field. The 
proposed solution relies on a distributed version of a filter, based on average 
consensus. Second, the forthcoming paper [10| develops a non parametric inter- 
polation algorithm: the field is defined on a compact subset of the plane, which 
is apportioned into regions assigned to each agent. Then, each node maintains 
a representation of its current region. Third, the papers [T2| [16] contain an 
articulated framework of cooperative exploration, in which a group of sensors 
travels in the plane, harvesting information about the field. 

With respect to these works, our paper differs on many respects. First, as we 
are interested in monitoring rather than in exploration, we choose to consider 
a static network, rather than a group of moving sensors. More important, we 
undertake one particular parametric estimation problem, which has not received 
any specific attention in the literature. This allows to obtain simple and neat 
results, and to focus on the core problem of fusing the local estimates obtained 
by the single nodes. 

1.1 Problem statement 

Let F : ^ M be a Gaussian function 

F{xi,X2) = Cie ^2 (1) 

where Ci and C2 are positive constants and m = (toi , TO2) is a point of the plane. 
Let us discuss the meaning of these parameters, considering the application of 
observing a diffusion process (fire, pollution). 

(i) The point m is the center of the distribution, that is the center of the fire 
or the location of the pollutant leak. 

(ii) The positive scalar Ci gives the amplitude of F, that is Ci = F{mi,m2). 

(iii) The positive scalar C2 gives information about the width of the event. If 
the distribution is the outcome of a diffusion process, C2 is proportional 
to the square of the time elapsed from the beginning of the event. 

Let there be a network of sensors, referred to as nodes. Each sensor is able 
to measure the value of F in its own location, and to communicate with its 
neighbors on the graph with no errors or delays. We want the network to be 
able to collectively estimate of the four parameters describing the Gaussian 
function. To this goal, we devise a two-phases algorithm. 



1) The sensors, communicating with their neighbors, construct local esti- 
mates of the parameters. 

2) The local estimates are fused into common estimates, using an iterative 
consensus algorithm. 

1.2 Statement of contributions and outline 

We design the two phases of distributed sensing and fusion separately. In Sec- 
tion [2] we discuss the geometric design of the sensor network, and we choose 
to work with an hexagonal tessellation. In Section [3] we consider how a single 
node of such network can build an estimate of the parameters of the Gaussian 
field, using local measurements, that is measurements taken by itself and its 
neighbors. Moreover, we study the the sensitivity to measurement errors of 
the estimates of the parameters, and the way the induced estimation error de- 
pends on the distance among sensors. Afterwards, in Section^ we consider the 
whole network, and we discuss how local estimates can be fused into a common 
estimate, through an iterative consensus algorithm. We show that fusion can 
effectively compensate errors in the measurements. Finally, Section [S] contains 
some concluding remarks. 

2 Network design: hexagonal tessellation 

In this section we discuss the properties for a network of sensors which are 
desirable to solve the problem stated in Section 11.11 We choose to work with 
networks whose nodes and links form (a portion of) a regular hexagonal tessel- 
lation of the plane. Indeed, an hexagonal tessellation enjoys many properties 
which make it the right candidate for our network. 

(i) Equispacing. In an hexagonal tessellation, all neighboring nodes are I 
far from each other, if we denote by / the length of the edge of each 
regular hexagon. If the communication among nodes is wireless, which 
is a very natural choice for a wide area network to be deployed in an 
environment, a common communication model is the disk graph: each 
node can communicate with nodes closer than a threshold R. Since in a 
regular hexagonal tessellation all nodes are I far apart from their neighbors, 
choosing I < R we satisfy the connectivity constraint. 

(ii) Modularity. Neglecting border effects, an hexagonal tessellation is vertex- 
transitive. This implies that the algorithmic design can easily be done 
off-line, and then be applied to any node. 

(iii) Coverage properties. Regular tessellations are optimal from the point of 
view of covering, since they are the Delaunay graph of a centroidal Voronoi 
partition [1]. 



In a distributed framework, every single node is supposed to obtain an estimate 
of the four parameters Ci, C2, mi, m2. We are going to prove in the sequel that, 
in order to do this, it needs at least four measurements, which can be its own 
measurement and three measurements from its neighbors. The graph of the 
hexagonal tessellation is 3-regular (every node has three neighbors), and then it 
allows the required information exchange. 

The above properties are satisfied not only by hexagonal tessellations, but 
also by the other regular tessellations of the plane, using triangles or squares. 
However, hexagons are preferable. Indeed, they require the minimum number of 
connections, three (instead of six and four, respectively): this can be a desirable 
property in terms of the communication cost. Moreover, hexagonal tessellations 
cover a larger area with the same number of nodes, and this is a useful feature, 
considering that the cost of putting in place a wireless network is proportional 
to the number of nodes. Indeed elementary formulas imply that, being N the 
number of nodes, the covered area is ^^PN for the hexagonal tessellation, while 
it is ^PN, and PN, for the triangular and square tessellations, respectively. 

Remark 2.1 The regularity properties of the network arc lost at the borders. 
For this reason, estimation is going to be performed at inner nodes only: such 
limitation is intuitively going to be less important for large networks. 



3 Estimation 

In this section, we show how an inner node of the hexagonal tessellation graph 
can estimate the four parameters of the Gausssian, using local information. 
Since all edges are of length /, up to isometrics we can assume the four sensors 
to be located in 




and labeled as in Figure [TJ 



Proposition 3.1 Let /ii, yU2, /is, ^4 be the measurements taken by the four 
sensors, labeled as in Figure [IJ Then the four parameters in ([T]) are given by 
the following formulas: 
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<i>(Ci,C2,mi,m2)= (^F{0,0),F{0,1),f[~^1,-^),f[^1,-!^^) (2) 
Then the problem consists in inverting $, that is solving the algebraic system 
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Finally, Ci = yUie . We have shown that, provided {(^1, /i2, /J.3, /i4) G 

I^>o • A*2M3/^4 < /^ili there is only one compatible Gaussian function, whose 
parameters are stated in the thesis. □ 



3.1 Estimation errors 



Let us assume that the measurements be affected by additive noise. This section 
is devoted to study the error committed computing the parameters from such 
noisy measurements. 

For i e {1, ... ,4}, let Aj be the error affecting measurement Hi and let A 
be the column vector having components A,;. If we recall the definition ([2]) of 
the function then the first order approximation of the error committed about 
the estimated parameters of F can be computed as (_D$)~^A. Let us assume 
that the measurement errors Ai be independent Gaussian random variables with 
E[Ai] = and common variance ]E[A|] = . Then the first order approxima- 
tions of the errors on the estimates of the parameters of F are Gaussian random 
variables: they have zero mean and their variance can be computed, after some 
algebra, as 
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Remark 3.2 (Center position and errors) Remark that all the above func- 
tions grow unbounded as |m| goes to infinity. This is consistent with the intu- 
ition that the sensors be more effective if close to the center of the Gaussian. 



Moreover, ^(Z; arctan ^) and S{l;\m\) grow unbounded also when |m| 0. 
This suggests the opportunity of having nodes enough far apart from each other. 

Hence, there is some non-trivial optimization issue about /. The problem of a 
local optimization of the inter-node distance is the topic of the next paragraph. 

3.2 Optimal spacing 

Since we have restricted ourselves to hexagonal regular networks, the design 
parameter we are left with for optimization is L It is clear that the choice of / 
is constrained by the area to be covered and the number of available sensors. 

Although we are aware of these potential constraints, in this section we keep 
them aside, and we we consider the problem of finding, for a fixed node, the 
value of / which minimizes each of the above errors variances. That is, we 
address the unconstrained optimization problem 

^opt = argmin5(0- (3) 

Whenever it is useful to specify which function we are referring to, we shall 
write /opt(C'i), Zopt(C'2), and so on. 

Remark 3.3 (Finding Zopt) Some facts on the optimization problem ([3]) are 
immediate. 

(i) The errors variances S{1; •) depend continuously on / and grow unbounded 
as I goes to or to oo. Hence there is at least one optimal I. Finding closed 
form expressions for these minima seems an intractable problem, but the 
minima can be evaluated numerically. 

(ii) The parameter Ci has no effect on any /opt- 

(iii) The parameter C2 can be set as C2 = 1 without loss of generality, since 
the mapping (m,C2) 1— > (Am, A^C2) induces a mapping as /opt — ^ A/optj 
for aU A e R>o. 

Hence m, the relative location of the center, is the most meaningful variable. 
Before showing some numerical results, we state a bound about S{l;C2)- its 
proof, which is based on direct computations, is postponed to the Appendix. 

Proposition 3.4 Consider the optimization problem Then 

VCh-lm] </opt(C2) < V2^+\m\. 

Namely, if m ~ 0, then /opt — 1.11691\/C2^. 

In Figures [2] and [31 we show some computational results about the minima. 
Some remarks are in order. 
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Figure 2: Numerical computation of Zopt(C2), depending on mi,m2, with 
C2 = 1. 



(i) The optimal Zopt depends on the parameters value, and namely on the 
value of m with respect to the central node, which is assumed to be located 
at the origin. As a consequence, the performed optimization is node- 
dependent, and can be of interest for off-line design when it is available 
some a priori information about the parameters. 

(ii) The functions S{1) have several local minima, and this results in disconti- 
nuities of the computed global minimum point ^opt- 

(iii) It appears that lim|,„|_>o2 Zopt = 0. However, we know from Remark 13.21 
that in such limit S{1) diverges, for any I. Hence the estimate becomes 
unreliable, and the optimization looses its significance. We shall see, in 
the next section, how a data fusion algorithm can effectively cope with 
this problem, focusing on sensors close to the Gaussian center, so that the 
optimization of I for far-away sensors becomes irrelevant. 



4 Fusion by consensus 

In this section we discuss how the estimates of the parameters, obtained lo- 
cally by the single nodes, can be fused at the network level into an agreement 
about the parameters of the Gaussian field. To the goal of fusion, we can take 
advantage of the wide literature [13] on distributed consensus problems. Such 
problems can be solved by iterative distributed algorithms as the one we de- 
scribe. 



Figure 3: Numerical computation of ^optdwl), depending on rni,m2, with 
C2 = 1. 



To work with more gcnerahty, let us consider a quantity of interest a; G M, 
which in the case of this paper wiU be one of the parameters of the Gaussian 
F. The algorithm evolves through discrete time steps t S Z>o. The scalar 
Xi{t) denotes the local estimate that, at time t, sensor i e {1, • • • ,N} has of 
the quantity x: the estimate obtained from the estimation phase, described in 
Section [3l will be Xi{0). Note that, in our apphcation, only the inner nodes 
of the network do possess initial estimates, and then take part to the fusion 
process. Thus N is the number of inner nodes in the network. The iteration of 
the algorithm consists in 

x{t+l) = P{t)xit), (4) 

where x{t) = [xi (i), . . . , XN{t)]'^ is the estimates vector, and P{t) is a stochastic 
matrix, depending on the network topology and on time. Namely, Pij(t) = if i 
and j are not neighbors. The algorithm above, under suitable assumptions [3], 
can be proved to converge to a consensus, in the sense that for alH e {1, . . . , N}, 
we have limt_>oo (t) = x. The value x is said to be the consensus value, and 
it is a convex combination of the initial values x(0). In some cases, the weights 
of this convex combination can be chosen by the designer. For instance, if P{t) 
is constantly equal to a doubly stochastic matrix P, then x is the arithmetic 
average x^ve = N''^ Y^iLi ^iW- 

Averaging is the first natural idea for fusion of the estimates obtained by 
each node. Indeed, elementary probability formulas imply that, given that each 
estimate Xi (0) be affected by an independent random additive error whose vari- 
ance is af, the average x^vc = J2iLi is affected by an error of variance 
7V~^ '^i- facts, if all the errors are equal, this yields reducing the error 



by a factor A''. However, this is not the case of the present paper, when the 
quality of the estimates varies among nodes. Indeed, it can very weh happen 
that, for some node j, the error variance be (t| < J2iLi '^f- We argue that 
computing the arithmetic average is not always the best choice. To cope with 
this problem wc need to use a better estimator than the average. Given the 
noisy estimates Xi{Q), such that E[(a;i(0) — x^] = a'f, the optimal estimator of 
a; in a mean squared error sense is 

If the variances for all /c e {!,..., N} were known a priori, then the problem 
could be solved as follows, combining two consensus algorithms. For every node 
fc e {1, . . . ,7V}, let afc(O) = and 6fc(0) = ^. Then, the time-dependent 

vectors a{t) and b{t) evolve following two standard average consensus algorithms. 
We obtain that, for all j € {1, ... , TV}, 

1 Y-JV xijO) 

aAt) N 2^1=1 
hm -r-rr = —, rr t- = x. 



Instead, the variances {ufc} arc not a priori known in the case of this paper, 
because they depend on the value x, as illustrated in Scction FS.ll For this reason, 
the above procedure can not be used, and we propose instead a suboptimal 
distributed algorithm. At each time step, every node i € {!,... , TV} updates 
and communicates both the variable Xi(t), and a variable Si{t) € M>o, which 
approximates the variance of the error that i is committing in estimating x with 
the value Xi{t). The algorithm we propose is as follows. 

Let us consider the graph Q of the connections among the inner nodes of 
the hexagonal network, whose nodes are labeled in {!,..., TV}. For all i G 
{1, . . . , N}, let Mi C {!,..., N} denote the set of the neighbors of i in the 
graph, including i itself. Let a;(0) be the vector of the local estimates, and s(0) 
the vector of the presumed variances of the errors affecting .t(0), as computed 
by each node i substituting Xi{0) into the formulas in Section [3. II 

At each time step t g Z>o, each node j communicates its states Xj{t) and 
Sj(t) to all its neighbors, and the states arc updated following 

x{t + 1) = P{t)x{t) (5) 
s{t + l)=M[t)s{t), 



where 



^ if J e A^i, 



otherwise; 



and 



if J eM 



,2 ' 







otherwise. 



The above update can be done in a distributed way, using information at 
the node and received from the immediate neighbors. Remark that the update 
matrices depend on the state of the system, and then on time. The convergence 
of the proposed algorithm is stated in the following result. 

Proposition 4.1 (Convergence) Provided the graph Q is connected, system ^ 
converges to consensus in both x{t) and s{t), that is, there exist scalars x* , s* 
such that 



where 1 is the N x 1 vector whose entries equal 1. 

Proof: Unlike time-independent consensus iterations [3], an eigenvalue 
analysis is not in general sufficient to conclude convergence for a time-varying 
system. Instead, we take advantage of the very general convergence results in 
[TT] . To this goal, we have to check that the linear iterations ([S]) satisfy the 
Strict Convexity Assumption in [iTl page 172]. We make the following remark: 
every M{t) is a stochastic matrix, and then Sj{t) G [min^ Si(0), max^ Si(0)], for 
all t > 0, and all j e {1, . . . , iV}. Since these intervals arc compact and do not 
contain zero, we argue that for all t > 0, the non-zero entries of P{t) and M{t) 
belong as well to a compact set of positive numbers. This is enough to satisfy 
the convexity assumption, and since Q is connected, we can apply [111 Th. 3] 
and conclude the proof. □ 
Note that, while the state space of each agent is two-dimensional, the object 
of interest is the consensus about the variable x{t) only. The consensus value x* 
is not in general equal to the optimal estimator x, but the simulations presented 
in Section [5] show that x* is much closer to the true value x, than the arithmetic 
average Xave- 

Remark 4.2 (Algorithm variations) The proposed algorithm is not the only 
heuristic solution: this remark informally illustrates two other possibilities. A 
very simple idea is the following one. Instead of running a consensus algorithm 
about s(t), the variables Si[t) can be computed at each time step from the cur- 
rent estimates Xi (t) , using the formulas in Section 13.11 A combination of this 
latter approach with a consensus scheme can also be pursued, as follows. At 
each time step t, each node i computes from Xi{t) an estimate Si{t] fc = 0), and 
then runs some (say k) consensus steps about s only. The obtained Si(t;k) is 
then used in a consensus step to compute Xi{t + 1). A convergence proof for 
these variations can be adapted from the one given above. 






Figure 4: Raw and fused estimates of the center of the Gaussian whose 
position (toi,TO2) is equal to (0,0), (0.5,0.5), (1, 1), and (1.5, 1.5), respectively. 
The label Wise consensus denotes the algorithm ([5]). 



5 Simulations 

Simulations arc presented to demonstrate the effectiveness of the proposed esti- 
mation and fusion algorithm. To provide an example, we focus on the estimation 
of the center (mi,m2) of the Gaussian field. We consider a network of twelve 
nodes, six of which are inner nodes forming an hexagon, and hence are able 
to compute estimates and participate to the fusion process. We assume that 
each measurement is corrupted by an additive random Gaussian variable whose 
variance is 1% of the maximum value of the signal to be estimated. Figure [4] 
shows the raw estimates obtained by each node, and the fused estimate of the 
center, compared with its true position. Note that, although the measurement 
error seems to be low, it is enough to perturb significantly the local estimates. 
Hence, the simulations setting is suitable to demonstrate the effectiveness of the 
proposed fusion algorithm. 

Simulations show that the estimation of the bi-dimensional parameter is 
possible with the proposed distributed method, and that the sensitivity to mea- 
surement errors is greatly reduced by fusion: the consensus algorithm proposed 
in ([5]) is much more accurate, and then preferable to simple averaging. As ex- 



pected, the sensor network is less effective when the Gaussian is eccentric with 
respect to the network. Nevertheless, for moderate eccentricity the fusion algo- 
rithm gives an effective compensation of the errors affecting the local estimates. 

6 Conclusion and open problems 

In this work we have posed and solved the problem of distributely estimating 
a Gaussian signal from noisy measurements taken by networked sensors. Our 
solution consists in a two-phases algorithm designed for a network which forms 
an hexagonal tessellation. In a first phase, sensors locally compute estimates of 
the parameters identifying the Gaussian signal, using their own measurement, 
and the measurements by their immediate neighbors. In the second phase, 
an heuristic consensus algorithm is run to obtain a common estimate. The 
consensus phase takes advantage of an estimation of the statistics of the errors 
affecting the local estimates. This feature allows to obtain encouraging results 
in terms of accuracy, also in presence of significant measurement noise, which 
yields large errors in the local estimates. Besides proving the convergence of 
the algorithm, we also studied the reliability of the local estimates, obtaining 
useful results for the network design. Our work is suitable of several further 
developments, that we briefly outline as follows. 

(i) The analysis given in the present paper can be extended, through similar 
computations, to the case of a more general function 

F{x) = Cie ^2 , 

provided F is a known positive definite matrix. This version of the problem 
has a clear interpretation in terms of the application to wild-fires mon- 
itoring, whenever there is an anisotropy in the fire spreading (e.g. due 
to prevalent winds). Beyond this extension, an approach similar to the 
one of the present paper is likely to be useful on more general parametric 
estimation problems. 

(ii) Assuming that the sensors be deployed to form an hexagonal network 
offers some advantages, but is also a restriction, compared to general sen- 
sor networks. This remark suggests some natural developments. First, a 
necessary step towards application is a sensitivity analysis with respect 
to errors in the deployment of sensors. Second, the application of the 
algorithm would benefit from the complementary design of distributed al- 
gorithms to achieve the hexagonal network by deployment of self-propelled 
sensors. Third, there can be cases in which the hexagonal configuration 
is difficult or impossible to achieve. Hence, it would be useful to extend 
our analysis to different, possibly non regular, graphs. Natural questions 
include: which graphs are suitable to solve the Gaussian parameters es- 
timation problem? Among those, which ones are preferable? After these 
basic questions are answered, a locational optimization problem could be 



defined and hopefully solved in a distributed fashion, to achieve an op- 
timal sensors deployment. Note that problems of distributed locational 
optimization have already been approached in the literature, for instance 
in HIT]. 

(iii) The literature on distributed sensor fusion has already devoted a consider- 
able attention to the design of algorithms to compute the average of local 
values. As we have shown, this can happen to be a non-optimal solution 
to a fusion problem, when local values have different (possibly unknown) 
reliability. Beyond any technical improvement of the proposed algorithm, 
we believe that further research should focus on the problem of an optimal 
distributed fusion of estimates having different "quality" . 
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A Proof of Proposition 13.41 

Up to a multiplicative constant, wc can compute 

d 4/ P-lm2-C2 
dl^^'^ ^ G2 ' 

P + ^Imi + l^ - C2 2!(!±:^!:!li±Z!2i 
+ 2 2 g2 

P - ^Imi + l^ - C2 . '('-v^"'i+"'2) 
H 2 2 ^ 

C2 

Hence, if to = 0, then the necessary condition reads — 1^ = 3, which 

can be numerically solved, giving I = 1.11691^6^. If instead m 7^ 0, then the 
minimum has to satisfy at least one of the following inequalities 

f -Im2-C2 > 

l^ + y^lm,+l^-C2 > (6) 



That is, 



lopt > min i ^ ( 7712 + \/ 4:0-2 + ml , , 



1712 


\/3toi 




2 






f + 


2 



/4c, + (!2J + ^f 



Such a minimum is always positive, and moreover 



2 
1 



4C, 



> 



m 



using the Cauchy-Schwartz inequahty and the fact that the function ?; i— > ?/ + 
a/ 46*2 + is increasing. 

Let us now look for an upper bound on /opt- To this goal, let us assume / to 
satisfy 



V -l\m\-C2 >0, 



(7) 



so that the inequalities ^ are all satisfied. Then, using > 1 + y, we argue 
that 



Sil) > l^ + P{\m\''-C2/2)~2Cl 



(8) 



Actually the right-hand-side of ([8]) is positive when 

f > i(C2/2- Imp + y'(C72/2-|mP)2 + 8C|). 



Since 



i (C2/2 - |m|2 + y^(C2/2-|7^|2)2+8C|) < v/|C2/2-|m|2|+y2C2 

< Iml + v/2C^ 



and, for every I > \m\ + ^26*2, Equation ([7]) is satisfied, we obtain the upper 

bound 

lopt < \m\ + \f2C2. 



